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Review  of  Objectives: 


The  research  undertaken  under  this  grant  has  been  focused  toward  performing  dynamic,  highly 
accurate  computational  analyses  of  ganged-inlet  configurations  (Project  1,  directed  by  Dr.  D.S. 
McRae)  and  hypersonic  engine  component  interactions  (Project  2,  directed  by  Dr.  J.R.  Edwards). 
The  specific  tasks  remain  the  same  as  indicated  in  the  original  proposal.  This  document 
summarizes  accomplishments  made  during  the  period  from  Aug.  31, 2000  to  September  1,  2001. 

II.  Status  of  Effort 

The  research  performed  during  this  period  has  produced  and  demonstrated  all  of  the  required 
codes,  algorithms  and  technology  needed  to  analyze  dynamically  scramjet  propulsion  modules, 
both  singly  and  ganged.  Cold-flow  unstart  of  a  generic  hydrocarbon-fuel  configuration  (tested  at 
NASA  Langley)  due  to  nozzle  backpressure  and  angle  of  attack  change  have  been  obtained  and 
animated  for  analysis.  An  inlet-isolator-combuster  hydrogen-fueled  configuration  tested  at  the 
National  Aerospace  Laboratory  in  Japan  has  been  simulated  dynamically  and  analysed  to 
determine  the  cause  of  flow  path  destabilization.  Finally,  an  algorithm  to  allow  dynamic 
adaptation  of  adjacent  blocks  has  been  developed  to  allow  resolution  of  all  of  the  physics  while 
maintaining  conservation  between  blocks.  The  codes  have  not  yet  been  integrated.  Follow-on 
funding  will  be  sought  to  do  so. 

II  Accomplishments: 

Project  1:  Ganged  Inlet  Analysis 

In  order  to  resolve  accurately  the  transients  that  occur  in  each  inlet/module,  a  dynamic  solution 
adaptive  mesh  algorithm  is  required  during  transient  calculations.  A  strongly-coupled  version  of 
this  algorithm  has  been  installed  in  a  perfect  gas  version  of  Edwards  highly  implicit  incomplete 
block  factorization  /  planar  Gauss-Seidel  algorithm.  This  algorithm  was  used  to  obtain  steady 
solutions  for  selected  conditions  and  geometiy  for  which  experimental  results  are  available.  Two 
different  unstart  modes  were  then  explored  with  this  code;  the  first  in  which  unstart  was  induced 
by  a  downstream  pressure  change  and  the  second  in  which  unstart  was  inducec  by  a  rapid  change 
in  inlet  angle-  of-  attack.  [1,2] 

However,  the  completion  of  these  examples  did  not  provide  all  of  the  capability  needed  to 
perform  well-  resolved  ganged  module  simulations.  In  order  to  illustrate  the  remaining  need,  it  is 
necessary  to  recognize  that  the  unstart  and  dynamic  mesh  adaptation  begins  in  one  module  prior 
to  any  fluid  communication  with  the  adjacent  module.  This  means  that  the  r-refinement  mesh' 
adjustment,  which  is  an  elliptical  process,  will  begin  to  move  mesh  nodes  toward  the  aft  section 
of  the  module/inlet  to  resolve  the  unstart  as  it  proceeds.  Since  each  inlet  grid  comprises  one  or 
more  individual  blocks,,  this  movement  will  require  mesh  nodes,  if  inter  block  grid  continuity  is 
maintained,  to  slide  onto  the  splitter  plate  between  the  modules,  presenting  a  difficult  logical 
problem.  An  obvious  solution  to  this  problem  is  to  allow  the  grid  to  become  discontinuous  as  the 
grid  adapts  in  one  block,  thereby  relieving  the  requirement  that  the  grid  in  the  adjacent  block  be 
slaved  to  the  adjacent  adaptation.  An  example  overlay  of  two  independently  adapted  block 
interface  grids  is  shown  in  Figure  (1).  The  usual  problem  encountered  in  this  case  is  that 
conservation  is  not  preserved  in  transferring  information  between  blocks.  A  technique  has  been 
developed  based  on  Ramshaws  algorithm  [3]  that  allows  independent  mesh  adaptation  in  each 
block  while  maintaining  fully  conservative  grid  block  interfaces.  [4,5] 


2 


1 

0.9 

0.8 

0.7 

0.6 

0.5 

0.4 

0.3 

0.2 

0.1 

0 


Figure  (1)  Example  of  discontinuous  block  interface  meshes  after  adaptation 


Many  multiblock  solvers  use  ghost  cells  (overlap  cells  taken  from  the  adjacent  block)  to  provide 
the  boundary  conditions  for  computations  within  each  block.  These  ghost  cells  are  chosen  or 
generated  (for  instance,  in  the  case  of  domain  decomposition)  so  that  they  can  provide  data  for 
the  computation  of  block  interface  cell  fluxes.  The  goal  is  to  have  the  flux  computation  proceed 
as  if  the  block  boundary  did  not  exist.  A  communication  algorithm  is  used  to  transfer  data  from 
the  interior  cells  of  the  adjacent  block  to  the  ghost  cells.  For  continuously  connected  blocks,  the 
communication  algorithm  involves  copying  dependent  variables  from  the  interior  cells  of  the 
adjacent  block  to  the  ghost  cells.  In  order  for  two  discontinuously  connected  adjacent  blocks  to 
communicate,  interior  data  of  the  adjacent  block  must  be  transferred  to  the  ghost  cells  using  a 
more  complex  algorithm.  For  two  3-D  blocks,  this  data  transfer  is  generally  a  3-D  problem  since 
the  pattern  of  overlap  between  the  ghost  cells  and  the  interior  cells  of  the  adjacent  block  vary  in 
all  three  dimensions.  However,  if  the  ghost  cells  are  constructed  appropriately,  this  data  transfer 
can  be  reduced  to  a  2-D  problem  in  computational  space  by  "tiling"  the  ghost  cells  from  the 
interface  cell  geometry.  An  explanation  of  the  technique  follows.  The  interface  of  two  3-D 
blocks  is  defined  by  two  surface  grids,  knoivn  as  the  block  boundary  grids,  which  are  defined  on 
the  boundaries  of  their  parent  blocks.  If  the  grid  blocks  are  continuous  across  the  interface,  then 
the  block  boundary  grids  are  identical.  If  the  grid  blocks  are  not  continuous  across  the  interface, 
then  a  four-sided  cell  face  belonging  to  the  boundaiy  grid  on  one  side  of  the  interface  is  divided 
into  N  polygons  by  the  boundary  grid  on  the  other  side  of  the  interface.  In  the  general  case,  a 
ghost  cell  would  be  cut  into  an  arbitrary  number  of  polyhedra  by  the  interior  cells  of  the  adjacent 
block,  but  if  ghost  cells  are  generated  properly,  then  the  ghost  cells  would  be  divided  into  N 
polyhedra,  just  as  their  faces  are  divided  into  N  polygons.  This  is  accomplished  by  locating  a 
boundaiy  point  in  the  (^,n)  transformed  computational  space  of  the  boundary  grid  on  the 
opposing  side  of  the  interface.  The  ghost  nodes  are  then  created  by  "tiling"  the  point,  which 
involves  using  the  same  ^  and  t)  coordinates  while  translating  by  ±  1  in  the  ^  direction.  The 
(^,  n,  0  coordinates  are  then  mapped  into  the  (x,y,z)  space  to  give  the  physical  location  of  the 
node.  Figure  (2)  depicts  the  results  of  this  process.[5] 
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Figure  (2)  Process  of  "Tiling"  to  Form  Ghost  Nodes 

Once  all  of  the  nodes  have  been  tiled,  the  cells  are  defined  by  the  nodes  using  either  a  cell- 
centered  or  node-centered  finite  volume  discretization.  If  ghost  cells  are  constructed  in  this  way, 
then  a  single  mapping  can  be  derived  which  transfers  dependent  variable  values  to  ghost  cells  of 
one  block  from  interior  cells  of  an  adjacent  block. .  As  previously  mentioned,  this  data  transfer 
is  simpler  if  the  ghost  nodes  are  tiled  in  order  to  generate  the  ghost  cells..  If  the  blocks 
are  labeled  one  and  two,  then  the  mapping  that  transfers  data  from  block  two  to  block  one  is 
/2_„, .  The  mapping  is  generally  sparse,  since  the  mapping  entry  linking  interior  cells  of 
block  two  to  ghost  cells  of  block  one  is  zero  except  when  the  cells  overlap. 

The  Ramshaw  algorithm  [3]  stems  from  the  general  equation  for  the  area  of  a  polygon  P: 

This  formula  allows  the  area  of  any  polygon  to  be  calculated  given  the  locations  of 
its  vertices.  Ramshaw  provides  a  systematic  version  for  accomplishing  this  process.  The  primary 
difficulty  is  that  of  determining  floating  point  line  coincidence  in  the  arbitrarily  overlaid  mesh.  To 
date,  all  of  the  problems  have  been  solved.  The  results  at  the  interface  for  a  advancing  front  that 
crosses  two  blocks  with  severely  discontinuous  interface  meshes  is  shown  in  Figure  (3).[4]  This 
algorithm  is  now  ready  for  inclusion  in  the  codes  described  in  the  following  paragraphs  and  those 
of  Project  2. 

The  technique  used  to  produce  the  unstarts  is  a  stand-alone  approach  for  structured  grid 
adaptation,  which  has  been  strongly  coupled  to  an  implicit  perfect-gas  solver  using  a  time- 
accurate,  dual  time-stepping  procedure.[6,7,8]  The  dual-time  stepping  approach  provides  a  self- 
consistent  means  of  correcting  the  flow  solution  to  correspond  with  the  newly-adapted  mesh.  The 
adaptation  algorithm  is  based  on  the  center-of-mass  concept  introduced  by  Benson  and  McRae 
[9]  and  is  applied  in  a  parametric  space,  essentially  defining  a  trilinear  interpolant  for  moving  grid 
points  in  physical  space  to  their  new  locations.  The  procedure  is  suitable  for  viscous  flows  on 
highly-stretched  meshes,  as  cell-aspect  ratio  and  cell  skewness  effects  are  incorporated  directly 
into  the  elliptic  equations  that  govern  the  solution  of  the  interpolation  variables. 
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Figure  (3)  Section  of  advancing  front  in  two  adjacent  blocks  with  discontinuous 

interface  grids. 


Defining  the  interpolation  variables  corresponding  to  the  computational  r\,  and  t;  directions  as 
p^,  p^,  and  p^ ,  elliptic  equations  of  the  following  form  are  solved  in  computational  space  to 

facilitate  adaptation: 
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In  this,  cOj^  are  weighting  functions  for  the  k*  interpolant.  The  major  difference  between  this 
approach  and  typical  center  of-mass  methods  lies  in  the  non-unity  coefficients  of  the  second  — 
derivative  terms.  These  are  functions  of  a  stretching  factor  ( |  VA:  1/|  | ,  for  example)  raised  to 

an  orthogonality  factor  (|  VA:  •  1  /(I  Vk  ||  |)  ).  These  coefficients  act  in  the  non-adapting 

directions  to  weight  the  point  movement  as  required  to  restore  acceptable  orthogonality. 
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Figure  4:  Adapted  mesh  within  scramjet  inlet  configuration 

Reference  [6]  presents  applications  of  this  concept  to  steady,  high-speed  inlet  flows,  while  [7,8] 
discuss  its  use  in  dynamic  simulations  of  unstart  of  a  Mach  3.5  scramjet  inlet.  A  snapshot  of  an 
adaptive  mesh  during  the  unstart  process  is  shown  in  Figure  4. 

Project  2:  Hypersonic  Engine  Component  Interactions 

Work  finalized  during  the  no-cost  extension  period  has  resulted  in  the  development  of  a  powerful 
procedure  for  conducting  unsteady  reactive-flow  calculations  on  massively-parallel 
machines.[10,l  1,12]  This  approach  combines  high  resolution  upwind  differencing  strategies  with 
a  dual-time  stepping  (or  subiteration)  procedure  for  recovering  second  order  temporal  accuracy. 

A  key  feature  of  the  approach  is  the  use  of  highly  implicit  incomplete  block  factorization  or 
planar  Gauss-Seidel  methods  to  alleviate  stability  restrictions  due  to  severe  grid  stretching.  This 
allows  the  use  of  physical  time  steps  much  larger  than  the  inviscid  stability  limit,  a  feature  that  is 
particularly  important  as  flowpath  responses  within  complete  engine  configurations  may  be  very 
slow,  compared  to  typical  characteristic  time  scales.  Computational  efficiency  is  maintained  by 
storing  the  factorization  of  system  Jacobian  matrix  in  core  memory  for  the  particular  block  (or 
group  of  blocks)  mapped  to  a  particular  processor.  After  initial  transients  have  been  purged,  the 
factorization  only  needs  to  be  re-evaluated  every  few  iterations,  significantly  reducing  the  overall 
expense.  Parallelization  of  the  solver  is  accomplished  through  standard  domain-decomposition 
strategies,  with  communication  between  processors  facilitated  by  MPI  routines.  Balakrishnan’s 
9  species  /  24  reaction  mechanism  [13]  is  currently  used  to  model  hydrogen  oxidation. 

Turbulent  effects  are  handled  either  by  the  Spalart-Allmaras  one-equation  turbulence  model  or 
Menter’s  hybrid  k-e  /  k-co  two-equation  turbulence  model.  A  more  recent  addition  is  a  capability 
for  blending  Menter’s  approach  for  near-wall  regions  with  a  large-eddy  simulation  technique  for 
ffee-shear  and  separated  flow  regions.[!4] 
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The  solver  has  been  validated  through  steady-state  simulations  of  the  3-D  shock  /  hydrogen  flame 
experiments  of  Driscoll  and  co-workers  [10],  among  other  cases.  Dynamic  simulations  of  the 
response  of  a  complete  scramjet  inlet-combustor  configuration  to  time-dependent  hydrogen  fuel 
injection  have  also  been  conducted  in  two  and  three  dimensions  [11,12].  These  simulations  have 
focused  on  resolving  the  mechanisms  leading  to  flowpath  destabilization  and  eventual  inlet 
unstart  due  to  excessive  heat  release.  Figure  5  presents  snapshots  of  temperature  contours  at 
different  spanwise  planes  of  a  scramjet  inlet-isolator-combustor  configuration,  tested  at  the 
National  Aerospace  Laboratory  (NAL),  Japan.  [15]  The  sequence  of  events  leading  to  flowpath 
destabilization  is  as  follows.  Hydrogen  fuel  is  entrained  into  pockets  of  separated  flow  located 
within  the  sidewall  boundary  layer.  The  fuel  then  mixes  and  ignites,  leading  to  flashback  (time 
t2).  The  upstream  propagation  and  growth  of  hot  regions  of  separated  flow  causes  significant 
blockage  of  the  inviscid  core  fluid,  raising  combustor  pressure  levels  and  eventually  resulting  in 
inlet  unstart  (time  ts).  This  mechanism  for  flowpath  destabilization  is  initiated  before  the 
combustor  becomes  thermally  choked  at  the  exit  plane. 

This  simulation  was  performed  on  a  mesh  containing  over  2  million  nodes,  with  15  equations 
solved  per  grid  point.  The  simulation  was  performed  on  40  processors  of  the  North  Carolina 
Supercomputing  Center’s  IBM  SP-2  and  required  about  4.5  days  of  run-time  to  complete. 


A:  centerine,  B:  quarter-plane,  C:  near  wall 


hydrogen  injection  through  stmt 
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Figure  5:  Temperature  contours  at  different  spanwise  planes  during  time-dependent 

scramjet  simulation 
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This  technology  is  currently  being  applied  to  the  GTX  rocket-based  combined-cycle  engine  [16] 
under  development  at  NASA  Glenn  Research  Center.  Research  will  be  focused  toward 
examining  the  dynamic  response  of  the  engine  flowfield  during  mode  transition. 

Work  to  be  performed  during  the  next  reporting  period 

A  final  report  will  be  delivered  to  AFOSR  as  per  grant  guidelines  and  reporting  requirements. 
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has  entered  the  Ph.D.  program  at  NCSU. 

Keith  S.  McDaniel  received  his  Master  of  Science  degree  in  May,  2001 .  He  is  currently 
employed  at  CFD  Research  Corporation,  Huntsville,  Al. 

V.  New  Discoveries 

None,  except  as  indicated  in  this  report. 
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